Analisi di Clustering sul Dataset Banknote
D’Atri, Lombardi, Maiello
04/04/2025
1. Introduzione
Il presente progetto si propone di applicare diverse tecniche di clustering al dataset banknote al fine di identificare gruppi omogenei presenti nei dati e confrontare l’efficacia dei vari metodi. L’analisi intende mettere in pratica le metodologie apprese nel modulo di Data Analysis & Statistical Learning, in particolare, verranno trattati i seguenti algoritmi di clustering: Clustering Gerarchico Agglomerativo e Divisivo, Clustering Partizionale (K-means e K-medoids) e Clustering Model-based (Mixture of Gaussian).
Tali algoritmo sono stati valutati su tre diverse tipologie di dataset: - Dataset Completo: utilizzo di tutte le variabili disponibili. - Dataset Ridotto: uso esclusivo delle variabili “Bottom” e “Diagonal”, individuate come le più discriminative. - Dataset Trasformato: applicazione della PCA per ridurre le variabili a componenti che sintetizzano l’informazione.
Per ciascuna configurazione, sono state valutate le performance dei metodi di clustering mediante metriche interne ed esterne. Il documento è organizzato nei seguenti capitoli:
- Introduzione: Presentazione del progetto, degli obiettivi e del contesto.
- Descrizione del Dataset e Analisi Esplorativa:
Esplorazione preliminare dei dati, analisi descrittiva e operazioni di
preprocessing, tra cui la costruzione del dataset trasformato mediante
la tecnica PCA.
- Applicazione delle Tecniche di Clustering: Implementazione dei metodi di clustering (gerarchico, partizionale e model-based).
- Validazione e Visualizzazione dei Risultati: Confronto dei risultati mediante metriche di validazione e visualizzazione dei cluster (con l’ausilio della PCA).
- Conclusioni e Sviluppi Futuri: Sintesi dei risultati ottenuti e prospettive per ulteriori analisi.
2. Descrizione e analisi del dataset
In questa sezione esamineremo il dataset banknote
per comprenderne la struttura, le caratteristiche e le relazioni tra le
variabili. Il dataset proviene dal pacchetto mclust e
contiene misure che possono essere utilizzate per distinguere tra
banconote autentiche e false (o per altri scopi diagnostici). Il dataset
contiene le seguenti variabili:
Status: lo stato della banconota (autentica o contraffatta)
Length: lunghezza della banconota (mm)
Left: larghezza del bordo sinistro (mm)
Right: larghezza del bordo destro (mm)
Bottom: larghezza del margine inferiore (mm)
Top: larghezza del margine superiore (mm)
Diagonal: lunghezza della diagonale (mm)
Di seguito carichiamo tutte le librerie necessarie per il progetto.
library("mclust")
library("reshape2")
library("ggplot2")
library("GGally")
library("dplyr")
library("clustertend")
library("FNN")
library("FactoMineR")
library("factoextra")
library("caret")
library("pheatmap")
library("clusterCrit")
library("NbClust")
library("fpc")
library("scatterplot3d")
library("corrplot")
library("cluster")
library("pROC")
library("gridExtra")Carichiamo il dataset.
data(banknote)Visualizziamo le prime righe del dataset e la sua struttura.
head(banknote)## Status Length Left Right Bottom Top Diagonal
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4str(banknote)## 'data.frame': 200 obs. of 7 variables:
## $ Status : Factor w/ 2 levels "counterfeit",..: 2 2 2 2 2 2 2 2 2 2 ...
## $ Length : num 215 215 215 215 215 ...
## $ Left : num 131 130 130 130 130 ...
## $ Right : num 131 130 130 130 130 ...
## $ Bottom : num 9 8.1 8.7 7.5 10.4 9 7.9 7.2 8.2 9.2 ...
## $ Top : num 9.7 9.5 9.6 10.4 7.7 10.1 9.6 10.7 11 10 ...
## $ Diagonal: num 141 142 142 142 142 ...Il dataset è composto da 200 osservazioni e 7 variabili, tutte di
tipo numerico ad eccezione di Status che è una variabile
binaria che può assumere uno tra i seguenti valori:
counterfeit e genuine.
Visualizziamo un riepilogo statistico delle variabili.
summary(banknote)## Status Length Left Right
## counterfeit:100 Min. :213.8 Min. :129.0 Min. :129.0
## genuine :100 1st Qu.:214.6 1st Qu.:129.9 1st Qu.:129.7
## Median :214.9 Median :130.2 Median :130.0
## Mean :214.9 Mean :130.1 Mean :130.0
## 3rd Qu.:215.1 3rd Qu.:130.4 3rd Qu.:130.2
## Max. :216.3 Max. :131.0 Max. :131.1
## Bottom Top Diagonal
## Min. : 7.200 Min. : 7.70 Min. :137.8
## 1st Qu.: 8.200 1st Qu.:10.10 1st Qu.:139.5
## Median : 9.100 Median :10.60 Median :140.4
## Mean : 9.418 Mean :10.65 Mean :140.5
## 3rd Qu.:10.600 3rd Qu.:11.20 3rd Qu.:141.5
## Max. :12.700 Max. :12.30 Max. :142.4Notiamo che la colonna Status è perfettamente
bilanciata.
Controlliamo se ci sono valori mancanti nelle colonne.
colSums(is.na(banknote))## Status Length Left Right Bottom Top Diagonal
## 0 0 0 0 0 0 0Come possiamo vedere nessuna colonna presenta valori mancanti. Per
avere una prima idea delle relazioni tra le variabili, creiamo una
matrice di scatterplot. Utilizziamo la colonna Status per
colorare i punti.
ggpairs(banknote[,-1],upper = list(continuous = "density", combo = "box_no_facet"),
lower = list(continuous = "points", combo = "dot_no_facet"),aes(colour=banknote$Status))
TODO: sintetizzare Lungo la diagonale sono visibili le distribuzioni
univariate delle variabili. Si nota una sovrapposizione parziale tra i
due gruppi di colore (Status), ma anche zone in cui uno dei due colori
prevale. Nei grafici di densità bivariata (metà diagonale superiore
della matrice) e negli scatter plot (metà diagonale inferiore della
matrice), si osservano alcune regioni in cui i punti di uno stesso
colore tendono a concentrarsi, suggerendo cluster parzialmente distinti.
Nonostante la sovrapposizione in diverse aree, si intravedono regioni in
cui un gruppo prevale. Ciò suggerisce che un algoritmo di
classificazione o di clustering potrebbe distinguere parzialmente le due
classi, anche se non in modo perfetto con una singola coppia di
variabili. È importante sottolineare che questo grafico mostra solo le
relazioni a coppie. Ciò significa che l’analisi si limita a considerare
due variabili alla volta, mentre una clusterizzazione netta potrebbe
emergere solo quando si considerano contemporaneamente più variabili. In
altre parole, anche se alcune coppie non evidenziano una separazione
chiara, combinando più dimensioni si potrebbe scoprire una struttura
clusterizzata più marcata. Il fatto che in alcune proiezioni (coppie di
variabili) si vedano zone di separazione o contorni distinti suggerisce
che esistono sottostrutture nei dati. Notiamo infine che le variabili
Bottom e Diagonal sembrerebbero separare in
modo migliore il dataset rispetto alle altre variabili.
D’ora in avanti, per ciascun algoritmo di clustering che proveremo, verranno effettuate tre diverse analisi:
Dataset completo: utilizzando tutte le variabili disponibili nel dataset.
Dataset ridotto: considerando esclusivamente le variabili
DiagonaleBottom.Dataset trasformato: Applicando l’Analisi delle Componenti Principali (PCA) per ridurre la dimensionalità.
Queste tre approcci permetteranno di valutare l’impatto delle diverse selezioni di variabili e della riduzione dimensionale sulle prestazioni e sui risultati degli algoritmi di clustering.
2.1. Valutazione della tendenza al clustering
Prima di applicare gli algoritmi di clustering, è fondamentale verificare se il dataset presenta una naturale propensione alla formazione di gruppi. Utilizziamo specifiche misure statistiche per determinare se i dati sono strutturati in cluster o se, al contrario, sono distribuiti in modo casuale. Questi metodi permettono di valutare se il clustering è effettivamente significativo e utile per l’analisi del dataset.
2.1.1. Hopkins statistics
L’Hopkins Statistic è una misura utilizzata per valutare la tendenza di un dataset alla clusterizzazione, confrontando la sua struttura con quella di un insieme di punti generati casualmente. Assume valori compresi tra 0 e 1, più il valore è vicino allo 0 più il dataset è strutturato in cluster.
banknote_num <- banknote[, -1]
banknote_num <- scale(banknote_num)
set.seed(987)
hopkins_stat <- hopkins(banknote_num, n = nrow(banknote_num) - 1)
random_df <- apply(banknote_num, 2, function(x){runif(length(x), min(x), max(x))})
random_df <- as.data.frame(random_df)
random_df <- scale(random_df)
set.seed(987)
hopkins_stat_rand <- hopkins(random_df, n = nrow(random_df)-1)
cat("Hopkins Statistics on Banknote dataset:", hopkins_stat$H, "\n")## Hopkins Statistics on Banknote dataset: 0.2642237cat("Hopkins Statistics on Random dataset:", hopkins_stat_rand$H, "\n")## Hopkins Statistics on Random dataset: 0.4926061Un valore di Hopkins pari a 0.2519264, essendo inferiore a 0.5, indica una tendenza alla clusterizzazione. In altre parole, le distanze tra i punti reali sono molto più ridotte rispetto a quelle dei punti generati casualmente, suggerendo la presenza di gruppi ben definiti nel dataset. Questo risultato supporta l’ipotesi che il dataset sia strutturato in cluster distinti.
2.1.2. Visual Assessment of cluster Tendency (VAT)
L’algoritmo VAT produce una matrice di dissimilarità ordinata, ottenuta riorganizzando le righe e le colonne della matrice delle distanze per evidenziare eventuali blocchi di punti “ravvicinati” (ossia potenziali cluster). Visualizziamo la matrice con l’intento di verificare se si formano lungo la diagonale secondaria del grafico, dei blocchi rossi che rappresentano regioni di osservazioni che condividono una forte similarità.
p1 <- fviz_dist(dist(banknote), show_labels = FALSE) + labs(title = "Banknote data")
p2 <- fviz_dist(dist(random_df), show_labels = FALSE) + labs(title = "Random data")
grid.arrange(p1, p2, ncol = 2, widths = c(6, 6))
Nei grafici VAT il rosso indica alta similarità (ossia bassa dissimilarità), mentre il blu indica bassa similarità (ossia alta dissimilarità). Nel grafico del dataset banknote, notiamo che lungo la diagonale sono presenti regioni rosse che indicano dati simili tra loro, suggerendo la presenza di cluster. Al contrario, il dataset randomico appare privo di struttura.
Possiamo concludere che il dataset Banknote ha una tendezza alla clusterizzazione, procediamo perciò con il calcolo della PCA e successivamente alla clusterizzazione.
2.2. Riduzione della dimensionalità con PCA
In questa sezione applichiamo la PCA sul dataset, considerando solo le colonne numeriche. Anche se le variabili sono tutte misurate in millimetri e quindi nella stessa scala, applichiamo comunque la standardizzazione per evitare che differenze nella dispersione (variabili con varianze maggiori) influenzino in modo sproporzionato i risultati della PCA.
banknote_num <- banknote[, -1]
head(banknote_num)## Length Left Right Bottom Top Diagonal
## 1 214.8 131.0 131.1 9.0 9.7 141.0
## 2 214.6 129.7 129.7 8.1 9.5 141.7
## 3 214.8 129.7 129.7 8.7 9.6 142.2
## 4 214.8 129.7 129.6 7.5 10.4 142.0
## 5 215.0 129.6 129.7 10.4 7.7 141.8
## 6 215.7 130.8 130.5 9.0 10.1 141.4pca_result <- PCA(banknote_num, scale = TRUE, graph = FALSE)Estraiamo gli autovalori e calcoliamo sia la varianza spiegata per ciascuna componente sia la varianza cumulativa. Queste informazioni ci serviranno per decidere il numero ottimale di componenti principali da mantenere.
eigenvalues <- get_eigenvalue(pca_result)[,1]
explained_variance <- eigenvalues / sum(eigenvalues)
cumulative_variance <- cumsum(explained_variance)
print(eigenvalues)## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6
## 2.9455582 1.2780838 0.8690326 0.4497687 0.2686769 0.1888799print(explained_variance)## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6
## 0.49092637 0.21301396 0.14483876 0.07496145 0.04477948 0.03147998print(cumulative_variance)## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6
## 0.4909264 0.7039403 0.8487791 0.9237405 0.9685200 1.0000000Abbiamo ora a disposizione gli autovalori e le informazioni sulla varianza, utili per valutare l’importanza di ciascuna componente principale.
Utilizziamo la Kaiser Rule, una tecnica utilizzata per determinare il numero di componenti principali da mantenere quando si esegue l’analisi delle componenti principali. La regola si basa sull’idea di conservare solo le componenti con una varianza spiegata maggiore di 1. Questo criterio si applica al valore proprio (eigenvalue) delle componenti principali.
kaiser_components <- sum(eigenvalues > 1)
cat("Numero di componenti selezionati con la Kaiser rule:", kaiser_components, "\n")## Numero di componenti selezionati con la Kaiser rule: 2Otteniamo 2 come numero di componenti principali da mantenere, ciò significa che le due componenti principali con valori propri superiori a 1 spiegano una quantità di varianza maggiore di quella che verrebbe spiegata da una singola variabile originale del dataset.
Un ulteriore criterio è quello di selezionare il numero minimo di componenti la cui somma di varianza spiegata (varianza comulativa), raggiunge almeno l’80% del totale. In questo modo garantiamo una buona rappresentazione dei dati.
variance_threshold <- min(which(cumulative_variance >= 0.80))
cat("Numero di componenti principali che spiegano l'80% della varianza comulativa:", variance_threshold, "\n")## Numero di componenti principali che spiegano l'80% della varianza comulativa: 3fviz_eig(pca_result, addlabels = TRUE, ylim = c(0, 50))
Nel grafico vediamo la varianza spiegata da ciascuna componente. Con le
prime tre componenti principali riusciamo a spiegare l’84.9% della
varianza cumulativa.
Vediamo infine un ultimo metodo per scegliere il numero di componenti principali: l’elbow method. Mediante un grafico scree plot visualizziamo gli autovalori per ciascuna componente e individuiamo “il gomito”, ovvero il punto nel grafico in cui la diminuzione degli autovalori rallenta. Tale punto può essere interpretato come un’indicazione del numero ottimale di componenti. Includiamo nel grafico anche due linee verticali per evidenziare i criteri della varianza cumulativa (rosso) e della Kaiser Rule (blu).
scree_plot <- data.frame(PC = 1:length(eigenvalues), Eigenvalue = eigenvalues)
ggplot(scree_plot, aes(x = PC, y = Eigenvalue)) +
geom_point() +
geom_line() +
geom_vline(xintercept = variance_threshold, linetype = "dashed", color = "red") +
geom_vline(xintercept = kaiser_components, linetype = "dotted", color = "blue") +
labs(title = "Scree Plot for PCA", x = "Principal Component", y = "Eigenvalue") +
theme_minimal()
Nel grafico il gomito sembra collocarsi intorno alla quarta
componente. I tre criteri danno perciò risultati differenti, per
mantenere un equilibrio tra interpretabilità e preservazione della
varianza, si potrebbe optare per 3 componenti.
Di seguito viene costruito il nuovo dataset in cui le variabili sono
sostituite con le prime 3 PCA, che verrà usato più avanti nelle
analisi.
pca_data <- as.data.frame(pca_result$ind$coord[, 1:3])
colnames(pca_data) <- c("PCA1", "PCA2", "PCA3")
banknote_pca <- cbind(pca_data, Status = banknote$Status)
banknote_pca <- banknote_pca %>% dplyr::select(Status, everything())
head(banknote_pca)## Status PCA1 PCA2 PCA3
## 1 genuine 1.7474012 1.65082829 1.4237611
## 2 genuine -2.2743177 -0.53879328 0.5326484
## 3 genuine -2.2774015 -0.10767707 0.7174149
## 4 genuine -2.2835546 -0.08765431 -0.6056336
## 5 genuine -2.6321283 0.03919590 3.1963847
## 6 genuine 0.7584073 3.08874513 0.7864804Per comprendere meglio il contributo delle variabili alle componenti principali, realizziamo una mappa di correlazione. Questo ci aiuta a capire come le variabili si relazionano alle componenti principali.
var <- get_pca_var(pca_result)
fviz_pca_var(pca_result, col.var = "contrib",
gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"))
Questo passaggio ci offre una visione immediata di come ogni variabile
contribuisce alle componenti principali e della loro interrelazione.
Osserviamo che le variaibili che appaiono con colori più accesi sono
Length che contribuisce molto a Dim.2, e le variabili
Diagonal, Left e Right con un
contributo importante a Dim.1 seppure opposto.
Approfondiamo l’analisi mostrando i contributi delle variabili attraverso diversi metodi grafici, in modo da identificare quali variabili influenzano maggiormente ciascuna componente.
corrplot(var$contrib, is.corr=FALSE)
L’analisi dei contributi evidenzia il peso relativo di ciascuna
variabile nelle componenti principali, utile per interpretare i
risultati della PCA. Si osserva che Length e
Top hanno contributi particolarmente elevati
rispettivamente a Dim.2 e Dim.3.
Approfondiamo come hanno contribuito le variabili per la prima componente principale.
fviz_contrib(pca_result, choice = "var", axes = 1)
Come era visibile già dal grafico precedente, le variabili che hanno
contribuito maggiormente alla prima componente principale sono
Diagonal, Right, Left.
Similmente al passaggio precedente, visualizziamo i contributi delle variabili per la seconda componente principale.
fviz_contrib(pca_result, choice = "var", axes = 2)
Notiamo che la variabile Length ha contribuito più del 60%
alla seconda componente principale.
3. Clustering gerarchico
3.1 Clustering Gerarchico Agglomerativo
Il Clustering Gerarchico Agglomerativo (HAC) è una tecnica di clustering che costruisce una gerarchia di cluster attraverso un approccio “bottom-up”. Inizia considerando ogni punto dati come un cluster separato e, successivamente, unisce iterativamente i cluster più simili fino a formare un unico cluster che racchiude tutti i dati. Questo metodo non richiede la specifica del numero di cluster a priori e produce una rappresentazione ad albero, nota come dendrogramma, che facilita la visualizzazione delle relazioni tra i cluster.
Per applicare il clustering gerarchico agglomerativo al nostro dataset, seguiamo questi passaggi:
Scaling: Standardizziamo i dati.
Calcolo della Matrice di Distanza: Determiniamo le distanze tra tutte le coppie di punti nel dataset. Una misura comune è la distanza euclidea.
Applicazione del clustering agglomerativo: Utilizziamo la funzione
hclust()per eseguire il clustering gerarchico e generare il dendrogramma.Taglio del Dendrogramma: Decidiamo il numero ottimale di cluster utilizzando varie tecniche, e tagliamo il dendrogramma al livello appropriato.
Valutazione: Valutiamo i risultati del clustering utilizzando varie metriche.
3.1.1. Dataset completo
Scaliamo il dataset e calcoliamo la relativa matrice di dissimilatità mediante distanza euclidea.
df <- scale(banknote_num)
distanze <- dist(df, method = "euclidean")
diss_matrix <- as.matrix(distanze)Applichiamo il clustering con metodo di linkage ward che
generalmente insieme all’average linkage, performa meglio.
hclust_clustering <- hclust(distanze, method = "ward.D2")Il dendrogramma mostra come le osservazioni si uniscono in cluster, con l’altezza dei rami che riflette la loro dissimilarità. Analizzando il dendrogramma, possiamo decidere il numero di cluster, tagliando l’albero a un’altezza che separa i gruppi in modo significativo.
3.1.1.1. Scelta del miglior numero di cluster k
fviz_dend(hclust_clustering, show_labels = FALSE, as.ggplot = TRUE)
Notiamo un ampio salto a circa 30 suggerisce una buona separazione in 2 cluster, mentre un altro salto meno marcato (intorno a 15) ne suggerisce 3. Per confermare il numero ottimale di cluster, si prosegue con metodi alternativi.
Osserviamo l’andamento del Total Within Sum of Squares (WSS) rispetto a k, cercando di individuare il punto in cui la curva inizia a “piegarsi” più marcatamente, cioè dove l’aggiunta di ulteriori cluster non porta a una riduzione significativa del WSS. Di seguito il grafico Elbow.
fviz_nbclust(df, FUN = hcut, method = "wss") +
geom_vline(xintercept = 2, linetype = 2) +
labs(subtitle = "Elbow method")
Tra k = 1 e k = 2, la diminuzione del WSS è molto forte, tra k = 2 e k =
3 c’è ancora un calo del WSS ma meno netto. Oltre k = 3, la pendenza
inizia a ridursi gradualmente. Dal grafico, sembra che il calo
principale avvenga per k = 2, mentre oltre si riduce più lentamente.
Analizziamo ora il grafico dell’Average Silhouette Width.
fviz_nbclust(df, FUN = hcut, method = "silhouette") +
labs(subtitle = "Silhouette method")
L’indice di silhouette risulta più alto in corrispondenza di k = 2. Questo suggerisce che la separazione tra i cluster e la coesione interna sia massima quando si dividono i dati in due gruppi. Da k = 3 l’indice diminuisce senza tornare ai valori di picco.
Osserviamo il grafico del Gap Statistic che confronta la dispersione all’interno dei cluster ottenuta dal clustering, con quella attesa in un dataset generato casualmente (quindi privo di gruppi).
set.seed(123)
fviz_nbclust(df, FUN = hcut, method = "gap_stat", nboot = 500, verbose = FALSE) +
labs(subtitle = "Gap statistic method")
Il grafico del Gap Statistic indica che quando il dataset viene diviso in 3 gruppi la struttura di clustering risulta migliore rispetto a quella attesa in un dataset casuale.
Utilizziamo infine il pacchetto NbClust che valuta il
numero ottimale di cluster utilizzando numerosi indici interni.
pdf(NULL)
nb <- NbClust(df, distance = "euclidean", min.nc = 2, max.nc = 10, method = "ward.D2")## *** : The Hubert index is a graphical method of determining the number of clusters.
## In the plot of Hubert index, we seek a significant knee that corresponds to a
## significant increase of the value of the measure i.e the significant peak in Hubert
## index second differences plot.
## ## *** : The D index is a graphical method of determining the number of clusters.
## In the plot of D index, we seek a significant knee (the significant peak in Dindex
## second differences plot) that corresponds to a significant increase of the value of
## the measure.
##
## *******************************************************************
## * Among all indices:
## * 7 proposed 2 as the best number of clusters
## * 10 proposed 3 as the best number of clusters
## * 1 proposed 4 as the best number of clusters
## * 1 proposed 7 as the best number of clusters
## * 1 proposed 8 as the best number of clusters
## * 3 proposed 10 as the best number of clusters
##
## ***** Conclusion *****
##
## * According to the majority rule, the best number of clusters is 3
##
##
## *******************************************************************dev.off()## png
## 2La decisione finale basata sulla majority rule suggerisce 3 cluster, indicando che, complessivamente, questa soluzione offre una migliore separazione e coesione interna rispetto alle altre soluzioni.
Complessivamente la metà dei metodi analizzati (gap statistic, nbclust) converge su una soluzione a 3 cluster, mentre l’altra metà (elbow, e silhouette) converge a 2 cluster. Pertanto non possiamo optare per una scelta di “maggioranza”, decidiamo di scegliere un numero di cluster pari a 2 poichè il salto più marcato nel dendrogramma indica una forte separazione naturale in 2 gruppi.
Procediamo con il taglio del dendogramma in modo da ottenere 2 cluster. Successivamente visualizziamo una tabella di frequenza che mostra il numero di osservazioni assegnate a ciascun cluster, permettendoci di verificare la distribuzione dei dati tra i 3 cluster ottenuti.
cluster_cut <- cutree(hclust_clustering, k = 2)
table(cluster_cut)## cluster_cut
## 1 2
## 102 98fviz_dend(hclust_clustering, k = 2,
cex = 0.5,
k_colors = c("#2E9FDF", "#E7B800"),
color_labels_by_k = TRUE,
rect = TRUE
)
3.1.1.2. Validazione
Valutiamo ora la qualità dei risultati del nostro algoritmo di clustering mediante quattro metriche: Average Silhouette Width, Dunn Index, Corrected Rand Index e Meila’s Variation of Information Index.
Iniziamo con la metrica della silhouette media (calcolata come la media di s(i) su tutte le osservazioni), che fornisce un’indicazione complessiva della qualità del clustering: si auspica un valore il più vicino ad 1 possibile.
sil <- silhouette(cluster_cut, distanze)
fviz_silhouette(sil, palette = "jco", ggtheme = theme_classic())## cluster size ave.sil.width
## 1 1 102 0.39
## 2 2 98 0.36
Un valore di 0.37 suggerisce che la struttura dei cluster non è
particolarmente forte: c’è una separazione moderata, ma c’è anche una
certa sovrapposizione tra i cluster.
Analizziamo ora il Dunn Index, che valuta la separazione tra i cluster rispetto alla loro compattezza, cioè quanto i cluster siano ben separati e internamente compatti. Se i dati contengono cluster compatti e ben separati, ci si aspetta che il diametro dei cluster sia piccolo e che la distanza tra i cluster sia grande. Pertanto l’indice di Dunn deve essere massimizzato.
dunn_index <- intCriteria(as.matrix(df), cluster_cut, "Dunn")$dunn
print(paste("Dunn Index:", dunn_index))## [1] "Dunn Index: 0.210791490416207"Un Dunn Index di 0.210 indica che i cluster hanno una separazione piuttosto bassa rispetto alla loro coesione interna, suggerendo possibili problemi di sovrapposizione.
Analizziamo il Corrected Rand Index che misura quanto bene il
clustering ottenuto corrisponde a una partizione di riferimento. Nel
nostro caso utilizziamo la colonna Status del dataset. Ed
infine il Meila’s Variation of Information (VI) Index che misura quanto
due cluster differiscono tra loro in termini di entropia. A differenza
del Corrected Rand Index, che valuta la similarità, il VI Index misura
la distanza tra due partizioni (va perciò minimizzato).
status <- as.numeric(banknote$Status)
clust_stats <- cluster.stats(d = dist(df), status, cluster_cut)
print(paste("Corrected Rand Index:", clust_stats$corrected.rand))## [1] "Corrected Rand Index: 0.960200080403878"print(paste("Meila's Index:", clust_stats$vi))## [1] "Meila's Index: 0.0982391266151992"Un Corrected Rand Index di 0.9602 è molto elevato, il che indica che il clustering ottenuto corrisponde molto bene alla partizione di riferimento. Inoltre, un Information Index di 0.098 è basso, suggerendo che la differenza informativa tra le partizioni è minima.
Concludendo ci potrebbe essere sovrapposizione tra i cluster, ma senza compromettere eccessivamente la qualità dell’assegnazione.
3.2.1.3. Ottimizzazione
Ottimizziamo il clustering esplorando diverse combinazioni di metodi di linkage e numero di cluster. Valutiamo ciascuna combinazione mediante le metriche viste prima. Successivamente, assegniamo un ranking per ciascuna metrica (invertendo il criterio per quelle da massimizzare) e sommiamo i punteggi. Infine, ordiniamo i risultati e selezioniamo la combinazione con il punteggio complessivo più basso come miglior clustering.
evaluate_clustering <- function(df, distanze, status, linkage_methods = c("ward.D2", "single", "complete", "average", "centroid"), min_k = 2, max_k = 10) {
results_all <- data.frame(Method = character(),
k = numeric(),
Silhouette = numeric(),
Dunn = numeric(),
Corrected_Rand = numeric(),
VI = numeric(),
stringsAsFactors = FALSE)
for(method in linkage_methods) {
hclust_obj <- hclust(distanze, method = method)
for(k in min_k:max_k) {
cluster_cut <- cutree(hclust_obj, k = k)
sil <- silhouette(cluster_cut, distanze)
avg_sil <- mean(sil[, 3])
dunn_index <- intCriteria(as.matrix(df), cluster_cut, "Dunn")$dunn
clust_stats <- cluster.stats(distanze, status, cluster_cut)
corrected_rand <- clust_stats$corrected.rand
vi_val <- clust_stats$vi
results_all <- rbind(results_all, data.frame(Method = method,
k = k,
Silhouette = avg_sil,
Dunn = dunn_index,
Corrected_Rand = corrected_rand,
VI = vi_val,
stringsAsFactors = FALSE))
}
}
return(results_all)
}
results_all <- evaluate_clustering(df, distanze, status)
best_silhouette <- results_all %>% group_by(Method) %>% filter(Silhouette == max(Silhouette))
print(best_silhouette)## # A tibble: 5 × 6
## # Groups: Method [5]
## Method k Silhouette Dunn Corrected_Rand VI
## <chr> <int> <dbl> <dbl> <dbl> <dbl>
## 1 ward.D2 2 0.375 0.211 0.960 0.0982
## 2 single 2 0.333 0.341 0 0.718
## 3 complete 2 0.346 0.125 0.591 0.566
## 4 average 2 0.333 0.341 0 0.718
## 5 centroid 2 0.333 0.341 0 0.718best_dunn <- results_all %>% group_by(Method) %>% filter(Dunn == max(Dunn))
print(best_dunn)## # A tibble: 5 × 6
## # Groups: Method [5]
## Method k Silhouette Dunn Corrected_Rand VI
## <chr> <int> <dbl> <dbl> <dbl> <dbl>
## 1 ward.D2 2 0.375 0.211 0.960 0.0982
## 2 single 2 0.333 0.341 0 0.718
## 3 complete 10 0.186 0.179 0.268 1.45
## 4 average 2 0.333 0.341 0 0.718
## 5 centroid 2 0.333 0.341 0 0.718best_corrected_rand <- results_all %>% group_by(Method) %>% filter(Corrected_Rand == max(Corrected_Rand))
print(best_corrected_rand)## # A tibble: 5 × 6
## # Groups: Method [5]
## Method k Silhouette Dunn Corrected_Rand VI
## <chr> <int> <dbl> <dbl> <dbl> <dbl>
## 1 ward.D2 2 0.375 0.211 0.960 0.0982
## 2 single 8 -0.0124 0.225 0.000320 0.896
## 3 complete 3 0.327 0.138 0.688 0.596
## 4 average 5 0.316 0.222 0.951 0.140
## 5 centroid 3 0.228 0.298 0.000101 0.742best_vi <- results_all %>% group_by(Method) %>% filter(VI == min(VI))
print(best_vi)## # A tibble: 5 × 6
## # Groups: Method [5]
## Method k Silhouette Dunn Corrected_Rand VI
## <chr> <int> <dbl> <dbl> <dbl> <dbl>
## 1 ward.D2 2 0.375 0.211 0.960 0.0982
## 2 single 2 0.333 0.341 0 0.718
## 3 complete 2 0.346 0.125 0.591 0.566
## 4 average 5 0.316 0.222 0.951 0.140
## 5 centroid 2 0.333 0.341 0 0.718results_ranked <- results_all %>%
mutate(
rank_sil = rank(-Silhouette, ties.method = "min"),
rank_dunn = rank(-Dunn, ties.method = "min"),
rank_corrected_rand = rank(-Corrected_Rand, ties.method = "min"),
rank_vi = rank(VI, ties.method = "min"),
total_rank = rank_sil + rank_dunn + rank_corrected_rand + rank_vi
)
best_overall <- results_ranked %>% arrange(total_rank) %>% head(1)
print(best_overall)## Method k Silhouette Dunn Corrected_Rand VI rank_sil rank_dunn
## 1 ward.D2 2 0.3748428 0.2107915 0.9602001 0.09823913 1 21
## rank_corrected_rand rank_vi total_rank
## 1 1 1 24Il risultato mostra che il clustering ottimale secondo queste metriche è Ward.D2 con 2 cluster, che è la combinazione che abbiamo testato precedentemente. La combinazione scelta restituisce i valori migliori possibili per le metriche Silhouette, Corrected Rand Index e VI (Rank 1), mentre risulta meno performante rispetto ad altri metodi, ma compensato dagli altri punteggi, per il Dunn Index (Rank 21).
Visualizziamo i cluster ottenuti.
fviz_cluster(list(data = df, cluster = cluster_cut), geom = "point", ellipse.type = "norm", palette="jco", ggtheme = theme_minimal())
3.1.2. Dataset trasformato
Vediamo come cambia il clustering se invece che considerare tutte le variabili del dataset, consideriamo solo le tre componenti principali trovate prima.
banknote_pca_num <- banknote_pca[,-1]
results_all_pca <- evaluate_clustering(banknote_pca_num, dist(banknote_pca_num), status)
results_ranked_pca <- results_all_pca %>%
mutate(
rank_sil = rank(-Silhouette, ties.method = "min"),
rank_dunn = rank(-Dunn, ties.method = "min"),
rank_corrected_rand = rank(-Corrected_Rand, ties.method = "min"),
rank_vi = rank(VI, ties.method = "min"),
total_rank = rank_sil + rank_dunn + rank_corrected_rand + rank_vi
)
best_overall_pca <- results_ranked_pca %>% arrange(total_rank) %>% head(1)
print(best_overall_pca)## Method k Silhouette Dunn Corrected_Rand VI rank_sil rank_dunn
## 1 average 2 0.431344 0.1056305 0.9020088 0.1997658 1 27
## rank_corrected_rand rank_vi total_rank
## 1 2 1 31TODO: sintetizzare Il risultato indica che, applicando il clustering agglomerativo con numero di cluster pari a 2 e metodo di linkage “average” al dataset ridotto a tre componenti principali, otteniamo la combinazione migliore in base ai 4 indici di valutazione, anche se i risultati presentano alcune peculiarità. Rispetto alle altre combinazioni testate, la silhouette e la VI hanno ottenuto i migliori punteggi (rank 1), il Corrected Rand è quasi il migliore (rank 2), invece l’indice Dunn, in questa combinazione, non performa bene rispetto ad altre configurazioni testate (rank 27). Rispetto al clustering agglomerativo applicato su tutte le variabili del dataset il clustering sul dataset ridotto mostra: - Un indice silhouette leggermente migliore (0.431) rispetto a quello su tutte le variabili (0.375), suggerendo che i punti sono, mediamente, meglio raggruppati all’interno dei cluster nel caso PCA. - Il clustering su tutte le variabili ha un indice Dunn migliore (0.211 vs 0.106), il che significa che la separazione minima tra cluster è maggiore e i cluster risultano più distinti. - Entrambe le soluzioni mostrano valori molto alti del Corrected Rand e bassi del VI, ma la soluzione con tutte le variabili è leggermente superiore (0.960 vs 0.902 e 0.098 vs 0.200), indicando una migliore aderenza ad una struttura di riferimento o una maggiore stabilità del clustering.
Mentre l’approccio con PCA (riduzione a 3 componenti) offre cluster con una buona compattezza (indice silhouette migliore), il clustering basato su tutte le variabili, utilizzando il metodo ward.D2, risulta complessivamente superiore per quanto riguarda la separazione tra cluster e la corrispondenza con una struttura di riferimento (migliori valori di Dunn, Corrected Rand e VI).
Vediamo graficamente come sono separati i due cluster.
final_pca <- eclust(banknote_pca_num, "hclust", k = 2, hc_metric = "euclidean", hc_method = "average", graph = FALSE)
clusters <- final_pca$cluster
scatterplot3d(banknote_pca_num,
color = clusters,
pch = 19,
main = "Clustering su PCA (PC1, PC2, PC3)",
xlab = "PC1", ylab = "PC2", zlab = "PC3")
legend("topright", legend = unique(clusters), col = unique(clusters), pch = 19, title = "Cluster")
3.1.3. Dataset ridotto
Vediamo infine come cambia il clustering considerando solo le
variabili Bottom e Diagonal che sembravano
suddividere meglio i dati in gruppi.
banknote_bottom_diagonal <- banknote %>% select(Bottom, Diagonal)
banknote_bd_scaled <- scale(banknote_bottom_diagonal)
results_all_bd <- evaluate_clustering(banknote_bd_scaled, dist(banknote_bd_scaled), status)
results_ranked_bd <- results_all_bd %>%
mutate(
rank_sil = rank(-Silhouette, ties.method = "min"),
rank_dunn = rank(-Dunn, ties.method = "min"),
rank_corrected_rand = rank(-Corrected_Rand, ties.method = "min"),
rank_vi = rank(VI, ties.method = "min"),
total_rank = rank_sil + rank_dunn + rank_corrected_rand + rank_vi
)
best_overall_bd <- results_ranked_bd %>% arrange(total_rank) %>% head(1)
print(best_overall_bd)## Method k Silhouette Dunn Corrected_Rand VI rank_sil rank_dunn
## 1 ward.D2 2 0.6245096 0.09556396 0.9602 0.1120031 1 11
## rank_corrected_rand rank_vi total_rank
## 1 2 2 16TODO: sintetizzare In questo caso la combinazione migliore è data da
metodo di linkage ward.D2 e numero cluster pari a 2.
Notiamo che in questo caso che il valore Silhouette 0.6245 è più alto
rispetto agli altri due approcci (0.431 per PCA e 0.375 usando tutte le
variabili), il che indica che i punti all’interno dei cluster sono ben
raggruppati e distinti l’uno dall’altro. Ciò supporta l’idea iniziale
che queste due variabili siano in grado di dividere bene il dataset. Il
Dunn index risulta più basso rispetto a quello ottenuto con tutte le
variabili (0.211) o con la PCA (0.106), il che potrebbe significare che
alcuni cluster risultano più vicini o con una certa dispersione interna.
Il Corrected Rand 0.9602 rimane molto elevato, segnalando una forte
coerenza con la struttura di riferimento attesa. L’indice VI è basso,
confermando una buona similarità o una discreta aderenza a una struttura
di riferimento.
Visualizziamo come il dataset è stato diviso in cluster.
final_bd <- eclust(banknote_bd_scaled, "hclust", k = 2, hc_metric = "euclidean", hc_method = "ward.D2", graph = FALSE)
fviz_cluster(final_bd, geom = "point", ellipse.type = "norm", palette="jco", ggtheme = theme_minimal())
4. Clustering Partizionale
TODO: inserire commento introduttivo
4.1. K-Means
TODO: inserire commento introduttivo
4.1.1. Dataset completo
TODO: inserire commento introduttivo
4.1.1.1. Scelta del miglior numero di cluster k
In questo capitolo i cluster vengono creati utilizzando l’algoritmo k-means, dopo aver scalato le variabili numeriche del dataset.
df <- scale(banknote[, -1])Successivamente si determina il numero ottimale di cluster da utilizzare nell’algoritmo, impiegando il metodo Elbow.
fviz_nbclust(df, kmeans, method = "wss") +
geom_vline(xintercept = 2, linetype = 2)+
labs(subtitle = "Elbow method")
Il metodo Silhouette.
fviz_nbclust(df, kmeans, method = "silhouette")+
labs(subtitle = "Silhouette method")
Ed infine la GAP analysis:
set.seed(123)
fviz_nbclust(df, kmeans, nstart = 25, method = "gap_stat", nboot = 50)+
labs(subtitle = "Gap statistic method")
Il primo metodo e il secondo indicano 2 come numero ottimale di cluster, il terzo ne suggerisce 3. Poiché il dataset è composto da 200 osservazioni, equamente suddivise tra banconote autentiche e contraffatte, si sceglie di impostare il numero di cluster a due e si procede con l’esecuzione dell’algoritmo.
4.1.1.2. Validazione
TODO: inserire questa parte, prendere spunto dal clustering gerarchico
4.1.2 Dataset ridotto
4.1.2.1. Scelta del miglior numero di cluster k
In questo capitolo i cluster vengono creati utilizzando l’algoritmo k-means, dopo aver scalato le variabili numeriche del dataset.
df_ridotto <- df[, c("Bottom", "Diagonal")]TODO: Successivamente si determina il numero ottimale di cluster da utilizzare nell’algoritmo, impiegando il metodo Elbow.
fviz_nbclust(df_ridotto, kmeans, method = "wss") +
geom_vline(xintercept = 4, linetype = 2)+
labs(subtitle = "Elbow method")
TODO: Il metodo Silhouette.
fviz_nbclust(df_ridotto, kmeans, method = "silhouette")+
labs(subtitle = "Silhouette method")
TODO: Ed infine la GAP analysis:
set.seed(123)
fviz_nbclust(df_ridotto, kmeans, nstart = 25, method = "gap_stat", nboot = 500)+
labs(subtitle = "Gap statistic method")
TODO: Il primo metodo indica quattro come numero ottimale di cluster, il secondo ne suggerisce due mentre la gap statistic tre. Poiché il dataset è composto da 200 osservazioni, equamente suddivise tra banconote autentiche e contraffatte, si sceglie di impostare il numero di cluster a due e si procede con l’esecuzione dell’algoritmo.
4.1.2.2. Validazione
TODO: da fare
4.1.3 Dataset trasformato
4.1.3.1. Scelta del miglior numero di cluster k
In questo capitolo i cluster vengono creati utilizzando l’algoritmo k-means, dopo aver scalato le variabili numeriche del dataset.
df_pca <- scale(banknote_pca[, -1])TODO: Successivamente si determina il numero ottimale di cluster da utilizzare nell’algoritmo, impiegando il metodo Elbow.
fviz_nbclust(df_pca, kmeans, method = "wss") +
geom_vline(xintercept = 4, linetype = 2)+
labs(subtitle = "Elbow method")
TODO: Il metodo Silhouette.
fviz_nbclust(df_pca, kmeans, method = "silhouette")+
labs(subtitle = "Silhouette method")
TODO: Ed infine la GAP analysis:
set.seed(123)
fviz_nbclust(df_pca, kmeans, nstart = 25, method = "gap_stat", nboot = 500, verbose=FALSE)+
labs(subtitle = "Gap statistic method")
TODO: Il primo metodo indica quattro come numero ottimale di cluster, il secondo ne suggerisce due mentre la gap statistic tre. Poiché il dataset è composto da 200 osservazioni, equamente suddivise tra banconote autentiche e contraffatte, si sceglie di impostare il numero di cluster a due e si procede con l’esecuzione dell’algoritmo.
4.1.3.2. Validazione
TODO: da fare
4.2 K-Medoids
TODO: da rivedere tutto
4.2.1 Dataset completo
4.2.1.1. Scelta del miglior numero di cluster k
Come fatto in precedenza per l’algoritmo k-means, in questo capitolo si cercherà di determinare il numero di cluster k utilizzando i metodi elbow. silhouette e della GAP analysis sfruttando la funzione pam.
Vediamo il metodo Elbow.
fviz_nbclust(df, pam, method = "wss") +
geom_vline(xintercept = 2, linetype = 2)+
labs(subtitle = "Elbow method")
Il metodo Silhouette.
fviz_nbclust(df, pam, method = "silhouette")+
labs(subtitle = "Silhouette method")
Ed infine la GAP analysis:
set.seed(123)
fviz_nbclust(df, pam, nstart = 25, method = "gap_stat", nboot = 50)+
labs(subtitle = "Gap statistic method")
Tutti metodi suggeriscono che il numero ottimale dei cluster è pari a 2.
library(cluster)
head(df, n = 3)## Length Left Right Bottom Top Diagonal
## [1,] -0.2549435 2.433346 2.829942 -0.2890067 -1.183765 0.4482473
## [2,] -0.7860757 -1.167507 -0.634788 -0.9120152 -1.432847 1.0557460
## [3,] -0.2549435 -1.167507 -0.634788 -0.4966762 -1.308306 1.4896737pam.res <- pam(df, 2)
print(pam.res)## Medoids:
## ID Length Left Right Bottom Top Diagonal
## [1,] 193 -0.5205096 0.4944248 0.6026155 0.9570103 0.5598130 -1.113892
## [2,] 47 -0.2549435 -0.6135300 -0.6347880 -0.7735689 -0.5610584 0.882175
## Clustering vector:
## [1] 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [38] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 2 2 2 2
## [75] 2 2 2 2 2 2 2 2 2 2 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 1 1 1 1 1 1 1 1 1 1
## [112] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [149] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [186] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## Objective function:
## build swap
## 1.842553 1.788867
##
## Available components:
## [1] "medoids" "id.med" "clustering" "objective" "isolation"
## [6] "clusinfo" "silinfo" "diss" "call" "data"dd <- cbind(banknote, cluster = pam.res$cluster)
head(dd, n = 8)## Status Length Left Right Bottom Top Diagonal cluster
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 2
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 2
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 2
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 2
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 2
## 7 genuine 215.5 129.5 129.7 7.9 9.6 141.6 2
## 8 genuine 214.5 129.6 129.2 7.2 10.7 141.7 2pam.res$medoids## Length Left Right Bottom Top Diagonal
## [1,] -0.5205096 0.4944248 0.6026155 0.9570103 0.5598130 -1.113892
## [2,] -0.2549435 -0.6135300 -0.6347880 -0.7735689 -0.5610584 0.882175head(pam.res$clustering)## [1] 1 2 2 2 2 2cl <- pam.res$clustering
pairs(df, gap=0, pch=cl, col=c("#00AFBB", "#FC4E07")[cl])
Come previsto, i cluster risultano molto simili.
fviz_cluster(pam.res,
palette = c("#00AFBB", "#FC4E07"),
ellipse.type = "t",
repel = TRUE,
ggtheme = theme_classic()
)
4.2.2 Dataset ridotto
4.2.2.1. Scelta del miglior numero di cluster k
Come fatto in precedenza per l’algoritmo k-means, in questo capitolo si cercherà di determinare il numero di cluster k utilizzando i metodi elbow. silhouette e della GAP analysis sfruttando la funzione pam.
Vediamo il metodo Elbow.
df_ridotto <- df[, c("Bottom", "Diagonal")]
fviz_nbclust(df_ridotto, pam, method = "wss") +
geom_vline(xintercept = 2, linetype = 2)+
labs(subtitle = "Elbow method")
Il metodo Silhouette.
fviz_nbclust(df_ridotto, pam, method = "silhouette")+
labs(subtitle = "Silhouette method")
Ed infine la GAP analysis:
set.seed(123)
fviz_nbclust(df_ridotto, pam, nstart = 25, method = "gap_stat", nboot = 50)+
labs(subtitle = "Gap statistic method")
df_ridotto <- df[, c("Bottom", "Diagonal")]
pam.res <- pam(df_ridotto, 2)
print(pam.res)## Medoids:
## ID Bottom Diagonal
## [1,] 47 -0.7735689 0.8821750
## [2,] 185 0.8877871 -0.8535357
## Clustering vector:
## [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [75] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 1 2 2 2 2 2 2 2 2
## [112] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [149] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [186] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## Objective function:
## build swap
## 0.8831314 0.6414212
##
## Available components:
## [1] "medoids" "id.med" "clustering" "objective" "isolation"
## [6] "clusinfo" "silinfo" "diss" "call" "data"dd <- cbind(banknote, cluster = pam.res$cluster)
head(dd, n = 8)## Status Length Left Right Bottom Top Diagonal cluster
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 1
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 1
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 1
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 1
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 1
## 7 genuine 215.5 129.5 129.7 7.9 9.6 141.6 1
## 8 genuine 214.5 129.6 129.2 7.2 10.7 141.7 1pam.res$medoids## Bottom Diagonal
## [1,] -0.7735689 0.8821750
## [2,] 0.8877871 -0.8535357head(pam.res$clustering)## [1] 1 1 1 1 1 1cl <- pam.res$clustering
pairs(df, gap=0, pch=cl, col=c("#00AFBB", "#FC4E07")[cl])
Come previsto, i cluster risultano molto simili.
fviz_cluster(pam.res,
palette = c("#00AFBB", "#FC4E07"),
ellipse.type = "t",
repel = TRUE,
ggtheme = theme_classic()
)
4.2.3 Dataset trasformato
4.2.3.1. Scelta del miglior numero di cluster k
Come fatto in precedenza per l’algoritmo k-means, in questo capitolo si cercherà di determinare il numero di cluster k utilizzando i metodi elbow. silhouette e della GAP analysis sfruttando la funzione pam.
Vediamo il metodo Elbow.
df_pca <- scale(banknote_pca[, -1])
fviz_nbclust(df_ridotto, pam, method = "wss") +
geom_vline(xintercept = 2, linetype = 2)+
labs(subtitle = "Elbow method")
Il metodo Silhouette.
fviz_nbclust(df_pca, pam, method = "silhouette")+
labs(subtitle = "Silhouette method")
Ed infine la GAP analysis:
set.seed(123)
fviz_nbclust(df_pca, pam, nstart = 25, method = "gap_stat", nboot = 50)+
labs(subtitle = "Gap statistic method")
pam.res <- pam(df_pca, 2)
print(pam.res)## Medoids:
## ID PCA1 PCA2 PCA3
## 198 198 0.8816755 -0.1289878 0.19323058
## 21 21 -1.1983358 0.0749874 0.05222386
## Clustering vector:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
## 1 2 2 2 2 1 2 2 2 1 1 2 2 2 2 2 2 2 2 2
## 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40
## 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
## 2 2 2 2 2 2 2 2 2 2 2 1 2 2 2 2 2 2 2 2
## 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80
## 2 2 2 2 2 1 2 2 2 1 2 2 2 2 2 2 2 2 2 2
## 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
## 2 2 2 2 1 2 2 2 2 2 2 2 2 2 2 2 1 2 2 2
## 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## Objective function:
## build swap
## 1.370309 1.329691
##
## Available components:
## [1] "medoids" "id.med" "clustering" "objective" "isolation"
## [6] "clusinfo" "silinfo" "diss" "call" "data"dd <- cbind(banknote, cluster = pam.res$cluster)
head(dd, n = 8)## Status Length Left Right Bottom Top Diagonal cluster
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 2
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 2
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 2
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 2
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 1
## 7 genuine 215.5 129.5 129.7 7.9 9.6 141.6 2
## 8 genuine 214.5 129.6 129.2 7.2 10.7 141.7 2pam.res$medoids## PCA1 PCA2 PCA3
## 198 0.8816755 -0.1289878 0.19323058
## 21 -1.1983358 0.0749874 0.05222386head(pam.res$clustering)## 1 2 3 4 5 6
## 1 2 2 2 2 1cl <- pam.res$clustering
pairs(df, gap=0, pch=cl, col=c("#00AFBB", "#FC4E07")[cl])
Come previsto, i cluster risultano molto simili.
fviz_cluster(pam.res,
palette = c("#00AFBB", "#FC4E07"),
ellipse.type = "t",
repel = TRUE,
ggtheme = theme_classic()
)
TODO: da rivedere tutto
4.2.1 Dataset completo
In questo capitolo i cluster vengono creati utilizzando l’algoritmo k-medoids con k impostato a 2.
Come in precedenza, si analizzano prima i grafici a dispersione e successivamente i due cluster.
head(df, n = 3)## Length Left Right Bottom Top Diagonal
## [1,] -0.2549435 2.433346 2.829942 -0.2890067 -1.183765 0.4482473
## [2,] -0.7860757 -1.167507 -0.634788 -0.9120152 -1.432847 1.0557460
## [3,] -0.2549435 -1.167507 -0.634788 -0.4966762 -1.308306 1.4896737pam.res <- pam(df, 2)
print(pam.res)## Medoids:
## ID Length Left Right Bottom Top Diagonal
## [1,] 193 -0.5205096 0.4944248 0.6026155 0.9570103 0.5598130 -1.113892
## [2,] 47 -0.2549435 -0.6135300 -0.6347880 -0.7735689 -0.5610584 0.882175
## Clustering vector:
## [1] 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [38] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 2 2 2 2
## [75] 2 2 2 2 2 2 2 2 2 2 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 1 1 1 1 1 1 1 1 1 1
## [112] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [149] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [186] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## Objective function:
## build swap
## 1.842553 1.788867
##
## Available components:
## [1] "medoids" "id.med" "clustering" "objective" "isolation"
## [6] "clusinfo" "silinfo" "diss" "call" "data"dd <- cbind(banknote, cluster = pam.res$cluster)
head(dd, n = 8)## Status Length Left Right Bottom Top Diagonal cluster
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 2
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 2
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 2
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 2
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 2
## 7 genuine 215.5 129.5 129.7 7.9 9.6 141.6 2
## 8 genuine 214.5 129.6 129.2 7.2 10.7 141.7 2pam.res$medoids## Length Left Right Bottom Top Diagonal
## [1,] -0.5205096 0.4944248 0.6026155 0.9570103 0.5598130 -1.113892
## [2,] -0.2549435 -0.6135300 -0.6347880 -0.7735689 -0.5610584 0.882175head(pam.res$clustering)## [1] 1 2 2 2 2 2cl <- pam.res$clustering
pairs(df, gap=0, pch=cl, col=c("#00AFBB", "#FC4E07")[cl])
Come previsto, i cluster risultano molto simili.
fviz_cluster(pam.res,
palette = c("#00AFBB", "#FC4E07"),
ellipse.type = "t",
repel = TRUE,
ggtheme = theme_classic()
)
4.2.2 Dataset ridotto
df_ridotto <- df[, c("Bottom", "Diagonal")]
pam.res <- pam(df_ridotto, 2)
print(pam.res)## Medoids:
## ID Bottom Diagonal
## [1,] 47 -0.7735689 0.8821750
## [2,] 185 0.8877871 -0.8535357
## Clustering vector:
## [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [75] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 1 2 2 2 2 2 2 2 2
## [112] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [149] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [186] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## Objective function:
## build swap
## 0.8831314 0.6414212
##
## Available components:
## [1] "medoids" "id.med" "clustering" "objective" "isolation"
## [6] "clusinfo" "silinfo" "diss" "call" "data"dd <- cbind(banknote, cluster = pam.res$cluster)
head(dd, n = 8)## Status Length Left Right Bottom Top Diagonal cluster
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 1
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 1
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 1
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 1
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 1
## 7 genuine 215.5 129.5 129.7 7.9 9.6 141.6 1
## 8 genuine 214.5 129.6 129.2 7.2 10.7 141.7 1pam.res$medoids## Bottom Diagonal
## [1,] -0.7735689 0.8821750
## [2,] 0.8877871 -0.8535357head(pam.res$clustering)## [1] 1 1 1 1 1 1cl <- pam.res$clustering
pairs(df, gap=0, pch=cl, col=c("#00AFBB", "#FC4E07")[cl])
Come previsto, i cluster risultano molto simili.
fviz_cluster(pam.res,
palette = c("#00AFBB", "#FC4E07"),
ellipse.type = "t",
repel = TRUE,
ggtheme = theme_classic()
)
4.2.3 Dataset trasformato
df_pca <- scale(banknote_pca[, -1])
pam.res <- pam(df_pca, 2)
print(pam.res)## Medoids:
## ID PCA1 PCA2 PCA3
## 198 198 0.8816755 -0.1289878 0.19323058
## 21 21 -1.1983358 0.0749874 0.05222386
## Clustering vector:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
## 1 2 2 2 2 1 2 2 2 1 1 2 2 2 2 2 2 2 2 2
## 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40
## 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
## 2 2 2 2 2 2 2 2 2 2 2 1 2 2 2 2 2 2 2 2
## 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80
## 2 2 2 2 2 1 2 2 2 1 2 2 2 2 2 2 2 2 2 2
## 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
## 2 2 2 2 1 2 2 2 2 2 2 2 2 2 2 2 1 2 2 2
## 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## Objective function:
## build swap
## 1.370309 1.329691
##
## Available components:
## [1] "medoids" "id.med" "clustering" "objective" "isolation"
## [6] "clusinfo" "silinfo" "diss" "call" "data"dd <- cbind(banknote, cluster = pam.res$cluster)
head(dd, n = 8)## Status Length Left Right Bottom Top Diagonal cluster
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 2
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 2
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 2
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 2
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 1
## 7 genuine 215.5 129.5 129.7 7.9 9.6 141.6 2
## 8 genuine 214.5 129.6 129.2 7.2 10.7 141.7 2pam.res$medoids## PCA1 PCA2 PCA3
## 198 0.8816755 -0.1289878 0.19323058
## 21 -1.1983358 0.0749874 0.05222386head(pam.res$clustering)## 1 2 3 4 5 6
## 1 2 2 2 2 1cl <- pam.res$clustering
pairs(df, gap=0, pch=cl, col=c("#00AFBB", "#FC4E07")[cl])
Come previsto, i cluster risultano molto simili.
fviz_cluster(pam.res,
palette = c("#00AFBB", "#FC4E07"),
ellipse.type = "t",
repel = TRUE,
ggtheme = theme_classic()
)
5. Clustering Model-based: Mixture of Gaussian
In questa sezione introduciamo l’applicazione del clustering basato
su misture di gaussiane, utilizzando la funzione Mclust del
pacchetto mclust. Questo approccio assume che i dati
siano generati da una combinazione di distribuzioni gaussiane e, tramite
l’algoritmo EM (Expectation-Maximization), stima i parametri di ciascuna
componente. Il vantaggio principale di questo metodo è la capacità di
modellare cluster con forme ellissoidali, permettendo di gestire anche
cluster con sovrapposizioni. Inoltre, Mclust effettua una
selezione automatica del numero ottimale di cluster e del modello di
covarianza (parsimonioso) in base al criterio BIC, fornendo
un’indicazione della bontà dell’adattamento e della complessità del
modello.
data(banknote)
dati <- banknote[, -1]
dati_scaled <- scale(dati)
gmm_result <- Mclust(dati_scaled)
summary(gmm_result)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust VEE (ellipsoidal, equal shape and orientation) model with 4 components:
##
## log-likelihood n df BIC ICL
## -1191.595 200 51 -2653.405 -2666.898
##
## Clustering table:
## 1 2 3 4
## 16 99 47 38plot(gmm_result, what = "BIC")
plot(gmm_result, what = "classification")
plot(gmm_result, what = "uncertainty")
plot(gmm_result, what = "density")
Nonostante sappiamo che il dataset banknote presenta due classi reali
(ad esempio, banconote fraudolente e non fraudolente), l’analisi con
Mclust ha individuato come miglior modello un numero di cluster pari a
3. Inoltre, il modello selezionato è di tipo VEE (ovvero ellissoidale,
con forma variabile, ma con volume ed orientamento uguali), in base al
valore più alto di BIC evidenziato nel plot. È interessante notare come
il BIC per 2 cluster risulti significativamente inferiore rispetto a
quello per un numero maggiore di cluster, probabilmente perché sono
state utilizzate tutte le variabili del dataset per costruire i
cluster.
Confrontiamo esplicitamente i valori di BIC e ICL per soluzioni con 2 e
3 cluster:
data(banknote)
dati <- banknote[, -1]
dati_scaled <- scale(dati)
gmm_result2 <- Mclust(dati_scaled, G=2)
summary(gmm_result2)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust EEV (ellipsoidal, equal volume and shape) model with 2 components:
##
## log-likelihood n df BIC ICL
## -1260.823 200 49 -2781.264 -2781.417
##
## Clustering table:
## 1 2
## 101 99gmm_result3 <- Mclust(dati_scaled, G=3)
summary(gmm_result3)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust VEE (ellipsoidal, equal shape and orientation) model with 3 components:
##
## log-likelihood n df BIC ICL
## -1212.796 200 43 -2653.42 -2653.74
##
## Clustering table:
## 1 2 3
## 16 99 85Dai sommari ottenuti si osserva che entrambi i modelli, per k = 2 e k
= 3, forniscono valori che suggeriscono una soluzione migliore per 3
cluster. I valori del BIC e, conseguentemente, anche quelli dell’ICL
(che include una penalizzazione per l’incertezza nelle assegnazioni),
indicano che la soluzione a 3 cluster è preferibile rispetto a quella a
2 cluster.
Valutiamo l’effetto della scelta delle variabili sul modello.
Utilizziamo invece solo le due variabili “Bottom” e “Diagonal”,
individuate in precedenza dagli scatterplot come le più
discriminanti:
data(banknote)
dati <- banknote[, c("Bottom", "Diagonal")]
dati_scaled <- scale(dati)
gmm_result <- Mclust(dati_scaled)
summary(gmm_result)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust EVE (ellipsoidal, equal volume and orientation) model with 3 components:
##
## log-likelihood n df BIC ICL
## -387.8563 200 13 -844.5908 -849.2887
##
## Clustering table:
## 1 2 3
## 100 16 84plot(gmm_result, what = "BIC")
plot(gmm_result, what = "classification")
plot(gmm_result, what = "uncertainty")
plot(gmm_result, what = "density")
Anche in questo caso si nota che il valore di BIC è notevolmente più
alto (migliore) rispetto al caso in cui sono state usate tutte le
variabili. L’algoritmo suggerisce ancora un numero ottimale di cluster
pari a 3, ma questa volta ha selezionato come modello parsimonioso il
modello EVE (ellissoidale, con volume uguale e orientamento variabile).
È importante osservare che in uno dei cluster sono presenti solamente 16
osservazioni, il che potrebbe indicare un gruppo di dati che,
probabilmente, sono state erroneamente assegnate a un cluster
separato.
Verifichiamo ora se i valori di ICL sono in accordo con quelli del BIC
confrontando esplicitamente le soluzioni per 2 e 3 cluster:
data(banknote)
dati <- banknote[, c("Bottom", "Diagonal")]
dati_scaled <- scale(dati)
gmm_result2 <- Mclust(dati_scaled, G=2)
summary(gmm_result2)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust VVE (ellipsoidal, equal orientation) model with 2 components:
##
## log-likelihood n df BIC ICL
## -417.7539 200 10 -888.491 -890.1796
##
## Clustering table:
## 1 2
## 99 101gmm_result3 <- Mclust(dati_scaled, G=3)
summary(gmm_result3)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust EVE (ellipsoidal, equal volume and orientation) model with 3 components:
##
## log-likelihood n df BIC ICL
## -387.8563 200 13 -844.5908 -849.2887
##
## Clustering table:
## 1 2 3
## 100 16 84In questo caso, i valori di BIC e ICL sono concordi: per un numero di
cluster pari a 3 si ottiene un BIC di circa -845 e un ICL di -849,
mentre per 2 cluster il BIC risulta attorno a -888 e l’ICL a -890.
Questo dimostra che anche l’ICL, che penalizza ulteriormente le
assegnazioni incerte, suggerisce che il modello ottimale è quello con 3
cluster.
Adesso, procederemo ad utilizzare solamente le variabili estratte
tramite la tecnica PCA.
dati <- banknote_pca[, -1]
dati_scaled <- scale(dati)
gmm_result <- Mclust(dati_scaled)
summary(gmm_result)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust EEV (ellipsoidal, equal volume and shape) model with 2 components:
##
## log-likelihood n df BIC ICL
## -764.5395 200 16 -1613.852 -1616.362
##
## Clustering table:
## 1 2
## 102 98plot(gmm_result, what = "BIC")
plot(gmm_result, what = "classification")
plot(gmm_result, what = "uncertainty")
plot(gmm_result, what = "density")
In questo caso si nota che il valore di BIC è più alto (migliore)
rispetto al caso in cui sono state usate tutte le variabili ma è più
basso del caso in cui sono state usate solamente le due variabili
selezionate. Questa volta l’algoritmo suggerisce un numero ottimale di
cluster pari a 2, selezionando come modello parsimonioso il modello EEV
(ellissoidale, con volume e orientamento uguale).
Verifichiamo ora se i valori di ICL sono in accordo con quelli del BIC
confrontando esplicitamente le soluzioni per 2 e 3 cluster:
gmm_result2 <- Mclust(dati_scaled, G=2)
summary(gmm_result2)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust EEV (ellipsoidal, equal volume and shape) model with 2 components:
##
## log-likelihood n df BIC ICL
## -764.5395 200 16 -1613.852 -1616.362
##
## Clustering table:
## 1 2
## 102 98gmm_result3 <- Mclust(dati_scaled, G=3)
summary(gmm_result3)## ----------------------------------------------------
## Gaussian finite mixture model fitted by EM algorithm
## ----------------------------------------------------
##
## Mclust VEE (ellipsoidal, equal shape and orientation) model with 3 components:
##
## log-likelihood n df BIC ICL
## -764.778 200 19 -1630.224 -1657.822
##
## Clustering table:
## 1 2 3
## 64 98 38Anche in questo caso, i valori di BIC e ICL sono concordi: per un numero di cluster pari a 3 si ottiene un BIC di circa -1630 e un ICL di -1657, mentre per 2 cluster il BIC risulta attorno a -1614 e l’ICL a -1616. Questo dimostra che anche l’ICL suggerisce che il modello ottimale è quello con 2 cluster.
6. Classificazione con gli algoritmi di cluster
In questa sezione per le tre tipologie di dataset (completo, ridotto
e trasformato), analizzeremo la qualità dei diversi metodi di
clustering, fissando il numero di cluster pari al numero di classi
reali, ovvero 2. Questo approccio ci permette di utilizzare tecniche di
validazione esterne, dato che disponiamo delle etichette reali delle
osservazioni (variabile Status).
Per valutare l’efficacia dei cluster identificati, utilizzeremo diverse
metriche, tra cui:
- Adjusted Rand Index (ARI) che misura il grado di concordanza tra il
clustering ottenuto e la classificazione reale. - Metriche di
valutazione proprie del machine learning supervisionato, come
accuratezza, specificità, recall, matrice di confusione e AUC, per
interpretare il clustering come un problema di classificazione
binaria.
6.1. Dataset completo
Nel seguente codice, applichiamo i diversi algoritmi di clustering al dataset originale con i migliori parametri calcolati in precedenza (fissando il numero di cluster a 2) e aggiungiamo al dataframe le colonne corrispondenti alle assegnazioni di ciascun metodo. Da precisare che in aggiunta agli algoritmi discussi in precedenza, confrontiamo anche l’algoritmo gerarchico divisivo che non è stato affrontato prima.
data("banknote")
data <- banknote[, -1]
data_scaled <- scale(data)
set.seed(123)
km_model <- kmeans(data_scaled, centers = 2, nstart = 25)
pam_model <- pam(data_scaled, k = 2, nstart = 25)
hc_model <- hclust(dist(data_scaled, method = "euclidean"), method = "ward.D2")
diana_model <- diana(data_scaled)
gmm_model <- Mclust(data_scaled, G = 2)
banknote$kmeans <- km_model$cluster
banknote$kmedoids <- pam_model$clustering
banknote$hierarchical_agglomerative <- cutree(hc_model, k = 2)
banknote$hierarchical_divisive <- cutree(as.hclust(diana_model), k = 2)
banknote$gmm <- gmm_model$classification
head(banknote)## Status Length Left Right Bottom Top Diagonal kmeans kmedoids
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 1 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 2 2
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 2 2
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 2 2
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 2 2
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 1 2
## hierarchical_agglomerative hierarchical_divisive gmm
## 1 1 1 2
## 2 2 2 2
## 3 2 2 2
## 4 2 2 2
## 5 2 2 2
## 6 2 2 2Il codice successivo implementa una funzione per riallineare le assegnazioni dei cluster alle classi reali in modo da rendere i confronti più accurati.
bestMapping <- function(true, pred) {
banknote <- banknote %>%
mutate(Status_numeric = recode(Status,
"counterfeit" = 1,
"genuine" = 2))
true <- banknote$Status_numeric
true <- as.character(true)
pred <- as.character(pred)
levels_true <- sort(unique(true))
cm1 <- confusionMatrix(factor(pred, levels = levels_true), factor(true, levels = levels_true))
acc1 <- cm1$overall["Accuracy"]
pred_swap <- ifelse(pred == levels_true[1], levels_true[2], levels_true[1])
cm2 <- confusionMatrix(factor(pred_swap, levels = levels_true), factor(true, levels = levels_true))
acc2 <- cm2$overall["Accuracy"]
if(acc2 > acc1) {
return(list(mapped = pred_swap, cm = cm2))
} else {
return(list(mapped = pred, cm = cm1))
}
}Nel blocco successivo, calcoliamo le metriche di valutazione per ciascun algoritmo di clustering. È importante notare che il valore di AUC ha un significato più chiaro per l’algoritmo GMM (Gaussian Mixture Model), poiché questo metodo fornisce direttamente le probabilità a posteriori di appartenenza ai cluster. Per gli altri algoritmi, che operano in modo hard (assegnando ciascun punto a un solo cluster senza una misura di incertezza), il calcolo dell’AUC è stato forzato attraverso:
- Per K-means e K-medoids, è stato definito uno score basato sulla
distanza dei punti dal centroide o dal medoide del rispettivo
cluster.
- Per gli algoritmi gerarchici, non è stato possibile adottare un criterio analogo per costruire uno score continuo, motivo per cui l’AUC non è stato calcolato (indicata come NA nella tabella).
# K-means
ari_kmeans <- adjustedRandIndex(banknote$Status, banknote$kmeans)
mapping_kmeans <- bestMapping(banknote$Status, banknote$kmeans)
cm_kmeans <- mapping_kmeans$cm
acc_kmeans <- cm_kmeans$overall["Accuracy"]
sens_kmeans <- cm_kmeans$byClass["Sensitivity"]
spec_kmeans <- cm_kmeans$byClass["Specificity"]
distances_km <- t(apply(data_scaled, 1, function(x) {
apply(km_model$centers, 1, function(cent) sqrt(sum((x - cent)^2)))
}))
score_kmeans <- distances_km[,2] / rowSums(distances_km)
roc_obj_kmeans <- roc(banknote$Status, score_kmeans)
auc_kmeans <- auc(roc_obj_kmeans)
# K-medoids
ari_kmedoids <- adjustedRandIndex(banknote$Status, banknote$kmedoids)
mapping_kmedoids <- bestMapping(banknote$Status, banknote$kmedoids)
cm_kmedoids <- mapping_kmedoids$cm
acc_kmedoids <- cm_kmedoids$overall["Accuracy"]
sens_kmedoids <- cm_kmedoids$byClass["Sensitivity"]
spec_kmedoids <- cm_kmedoids$byClass["Specificity"]
medoids <- as.matrix(pam_model$medoids)
distances_kmed <- t(apply(data_scaled, 1, function(x) {
apply(medoids, 1, function(med) sqrt(sum((x - med)^2)))
}))
score_kmedoids <- distances_kmed[,2] / rowSums(distances_kmed)
roc_obj_kmedoids <- roc(banknote$Status, score_kmedoids)
auc_kmedoids <- auc(roc_obj_kmedoids)
# Gerarchico Agglomerativo
ari_hierAgg <- adjustedRandIndex(banknote$Status, banknote$hierarchical_agglomerative)
mapping_hierAgg <- bestMapping(banknote$Status, banknote$hierarchical_agglomerative)
cm_hierAgg <- mapping_hierAgg$cm
acc_hierAgg <- cm_hierAgg$overall["Accuracy"]
sens_hierAgg <- cm_hierAgg$byClass["Sensitivity"]
spec_hierAgg <- cm_hierAgg$byClass["Specificity"]
auc_hierAgg <- NA
# Gerarchico Divisivo
ari_hierDiv <- adjustedRandIndex(banknote$Status, banknote$hierarchical_divisive)
mapping_hierDiv <- bestMapping(banknote$Status, banknote$hierarchical_divisive)
cm_hierDiv <- mapping_hierDiv$cm
acc_hierDiv <- cm_hierDiv$overall["Accuracy"]
sens_hierDiv <- cm_hierDiv$byClass["Sensitivity"]
spec_hierDiv <- cm_hierDiv$byClass["Specificity"]
auc_hierDiv <- NA
# GMM
ari_gmm <- adjustedRandIndex(banknote$Status, banknote$gmm)
mapping_gmm <- bestMapping(banknote$Status, banknote$gmm)
cm_gmm <- mapping_gmm$cm
acc_gmm <- cm_gmm$overall["Accuracy"]
sens_gmm <- cm_gmm$byClass["Sensitivity"]
spec_gmm <- cm_gmm$byClass["Specificity"]
score_gmm <- gmm_model$z[,2]
roc_obj_gmm <- roc(banknote$Status, score_gmm)
auc_gmm <- auc(roc_obj_gmm)Infine, raccogliamo i risultati in un dataframe riassuntivo.
results_df <- data.frame(
Metodo = c("K-means", "K-medoids", "Gerarchico Agglomerativo", "Gerarchico Divisivo", "GMM"),
ARI = c(ari_kmeans, ari_kmedoids, ari_hierAgg, ari_hierDiv, ari_gmm),
Accuracy = c(acc_kmeans, acc_kmedoids, acc_hierAgg, acc_hierDiv, acc_gmm),
Sensitivity = c(sens_kmeans, sens_kmedoids, sens_hierAgg, sens_hierDiv, sens_gmm),
Specificity = c(spec_kmeans, spec_kmedoids, spec_hierAgg, spec_hierDiv, spec_gmm),
AUC = c(auc_kmeans, auc_kmedoids, auc_hierAgg, auc_hierDiv, auc_gmm)
)
print(results_df)## Metodo ARI Accuracy Sensitivity Specificity AUC
## 1 K-means 0.8456292 0.960 1.00 0.92 0.9973
## 2 K-medoids 0.9406018 0.985 1.00 0.97 0.9991
## 3 Gerarchico Agglomerativo 0.9602001 0.990 1.00 0.98 NA
## 4 Gerarchico Divisivo 0.9406015 0.985 0.99 0.98 NA
## 5 GMM 0.9799995 0.995 1.00 0.99 0.9999I risultati mostrano che il metodo GMM ha la migliore performance complessiva, con il più alto valore di ARI (0.98) e un’accuratezza del 99.5%. Il clustering gerarchico agglomerativo ha anch’esso buone prestazioni, mentre K-means e K-medoids mostrano una buona affidabilità ma con AUC leggermente inferiori. Analizzando la sensibilità e la specificità, si nota che la sensibilità – ovvero la capacità di identificare correttamente le banconote contraffatte – è generalmente molto alta per tutti gli algoritmi, indicando un’elevata capacità di individuare i falsi. Tuttavia, la specificità, ossia la capacità di riconoscere correttamente le banconote genuine, varia maggiormente tra i metodi. In particolare, K-means mostra una specificità inferiore (0.92) rispetto a K-medoids (0.97) e ai metodi gerarchici (0.98), suggerendo una maggiore tendenza a classificare erroneamente alcune banconote genuine come contraffatte; questo errore è meno grave di classificare una banconota contraffatta come vera.
6.2. Dataset ridotto
Nel seguente codice, applichiamo i diversi algoritmi di clustering al
dataset ridotto (con le sole variabili Diagonal e
Bottom) con i migliori parametri calcolati in precedenza
(fissando il numero di cluster a 2) e aggiungiamo al dataframe le
colonne corrispondenti alle assegnazioni di ciascun metodo.
data("banknote")
data <- banknote[, c("Bottom", "Diagonal")]
data_scaled <- scale(data)
set.seed(123)
km_model <- kmeans(data_scaled, centers = 2, nstart = 25)
pam_model <- pam(data_scaled, k = 2, nstart = 25)
hc_model <- hclust(dist(data_scaled, method = "euclidean"), method = "ward.D2")
diana_model <- diana(data_scaled)
gmm_model <- Mclust(data_scaled, G = 2)
banknote$kmeans <- km_model$cluster
banknote$kmedoids <- pam_model$clustering
banknote$hierarchical_agglomerative <- cutree(hc_model, k = 2)
banknote$hierarchical_divisive <- cutree(as.hclust(diana_model), k = 2)
banknote$gmm <- gmm_model$classification
head(banknote)## Status Length Left Right Bottom Top Diagonal kmeans kmedoids
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 2 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 2 1
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 2 1
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 2 1
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 2 1
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 2 1
## hierarchical_agglomerative hierarchical_divisive gmm
## 1 1 1 1
## 2 1 1 1
## 3 1 1 1
## 4 1 1 1
## 5 1 1 1
## 6 1 1 1Nel blocco successivo, calcoliamo le metriche di valutazione come nel caso precedente.
# K-means
ari_kmeans <- adjustedRandIndex(banknote$Status, banknote$kmeans)
mapping_kmeans <- bestMapping(banknote$Status, banknote$kmeans)
cm_kmeans <- mapping_kmeans$cm
acc_kmeans <- cm_kmeans$overall["Accuracy"]
sens_kmeans <- cm_kmeans$byClass["Sensitivity"]
spec_kmeans <- cm_kmeans$byClass["Specificity"]
distances_km <- t(apply(data_scaled, 1, function(x) {
apply(km_model$centers, 1, function(cent) sqrt(sum((x - cent)^2)))
}))
score_kmeans <- distances_km[,2] / rowSums(distances_km)
roc_obj_kmeans <- roc(banknote$Status, score_kmeans)
auc_kmeans <- auc(roc_obj_kmeans)
# K-medoids
ari_kmedoids <- adjustedRandIndex(banknote$Status, banknote$kmedoids)
mapping_kmedoids <- bestMapping(banknote$Status, banknote$kmedoids)
cm_kmedoids <- mapping_kmedoids$cm
acc_kmedoids <- cm_kmedoids$overall["Accuracy"]
sens_kmedoids <- cm_kmedoids$byClass["Sensitivity"]
spec_kmedoids <- cm_kmedoids$byClass["Specificity"]
medoids <- as.matrix(pam_model$medoids)
distances_kmed <- t(apply(data_scaled, 1, function(x) {
apply(medoids, 1, function(med) sqrt(sum((x - med)^2)))
}))
score_kmedoids <- distances_kmed[,2] / rowSums(distances_kmed)
roc_obj_kmedoids <- roc(banknote$Status, score_kmedoids)
auc_kmedoids <- auc(roc_obj_kmedoids)
# Gerarchico Agglomerativo
ari_hierAgg <- adjustedRandIndex(banknote$Status, banknote$hierarchical_agglomerative)
mapping_hierAgg <- bestMapping(banknote$Status, banknote$hierarchical_agglomerative)
cm_hierAgg <- mapping_hierAgg$cm
acc_hierAgg <- cm_hierAgg$overall["Accuracy"]
sens_hierAgg <- cm_hierAgg$byClass["Sensitivity"]
spec_hierAgg <- cm_hierAgg$byClass["Specificity"]
auc_hierAgg <- NA
# Gerarchico Divisivo
ari_hierDiv <- adjustedRandIndex(banknote$Status, banknote$hierarchical_divisive)
mapping_hierDiv <- bestMapping(banknote$Status, banknote$hierarchical_divisive)
cm_hierDiv <- mapping_hierDiv$cm
acc_hierDiv <- cm_hierDiv$overall["Accuracy"]
sens_hierDiv <- cm_hierDiv$byClass["Sensitivity"]
spec_hierDiv <- cm_hierDiv$byClass["Specificity"]
auc_hierDiv <- NA
# GMM
ari_gmm <- adjustedRandIndex(banknote$Status, banknote$gmm)
mapping_gmm <- bestMapping(banknote$Status, banknote$gmm)
cm_gmm <- mapping_gmm$cm
acc_gmm <- cm_gmm$overall["Accuracy"]
sens_gmm <- cm_gmm$byClass["Sensitivity"]
spec_gmm <- cm_gmm$byClass["Specificity"]
score_gmm <- gmm_model$z[,2]
roc_obj_gmm <- roc(banknote$Status, score_gmm)
auc_gmm <- auc(roc_obj_gmm)Infine, raccogliamo i risultati in un dataframe riassuntivo.
results_df <- data.frame(
Metodo = c("K-means", "K-medoids", "Gerarchico Agglomerativo", "Gerarchico Divisivo", "GMM"),
ARI = c(ari_kmeans, ari_kmedoids, ari_hierAgg, ari_hierDiv, ari_gmm),
Accuracy = c(acc_kmeans, acc_kmedoids, acc_hierAgg, acc_hierDiv, acc_gmm),
Sensitivity = c(sens_kmeans, sens_kmedoids, sens_hierAgg, sens_hierDiv, sens_gmm),
Specificity = c(spec_kmeans, spec_kmedoids, spec_hierAgg, spec_hierDiv, spec_gmm),
AUC = c(auc_kmeans, auc_kmedoids, auc_hierAgg, auc_hierDiv, auc_gmm)
)
print(results_df)## Metodo ARI Accuracy Sensitivity Specificity AUC
## 1 K-means 0.9799995 0.995 0.99 1.00 0.9999
## 2 K-medoids 0.9799995 0.995 0.99 1.00 0.9999
## 3 Gerarchico Agglomerativo 0.9602000 0.990 0.99 0.99 NA
## 4 Gerarchico Divisivo 0.9406015 0.985 0.98 0.99 NA
## 5 GMM 0.9799995 0.995 1.00 0.99 0.9995Questa configurazione semplifica il modello senza perdere
accuratezza, dimostrando che Bottom e Diagonal
sono le variabili chiave per distinguere le classi:
- K-means, K-medoids e GMM ottengono ARI 0.98, accuracy 99.5% e AUC ~1,
confermando che queste due variabili sono altamente
discriminative.
- GMM raggiunge sensibilità 100%, identificando tutte le banconote
contraffatte, mentre K-means e K-medoids garantiscono specificità 100%,
evitando falsi positivi.
6.3. Dataset ridotto
Nel seguente codice, applichiamo i diversi algoritmi di clustering al dataset trasformato (con le variabili trasformate mediante tecnica PCA) con i migliori parametri calcolati in precedenza (fissando il numero di cluster a 2) e aggiungiamo al dataframe le colonne corrispondenti alle assegnazioni di ciascun metodo.
data <- banknote_pca[, -1]
data_scaled <- scale(data)
set.seed(123)
km_model <- kmeans(data_scaled, centers = 2, nstart = 25)
pam_model <- pam(data_scaled, k = 2, nstart = 25)
hc_model <- hclust(dist(data_scaled, method = "euclidean"), method = "ward.D2")
diana_model <- diana(data_scaled)
gmm_model <- Mclust(data_scaled, G = 2)
banknote$kmeans <- km_model$cluster
banknote$kmedoids <- pam_model$clustering
banknote$hierarchical_agglomerative <- cutree(hc_model, k = 2)
banknote$hierarchical_divisive <- cutree(as.hclust(diana_model), k = 2)
banknote$gmm <- gmm_model$classification
head(banknote)## Status Length Left Right Bottom Top Diagonal kmeans kmedoids
## 1 genuine 214.8 131.0 131.1 9.0 9.7 141.0 2 1
## 2 genuine 214.6 129.7 129.7 8.1 9.5 141.7 1 2
## 3 genuine 214.8 129.7 129.7 8.7 9.6 142.2 1 2
## 4 genuine 214.8 129.7 129.6 7.5 10.4 142.0 1 2
## 5 genuine 215.0 129.6 129.7 10.4 7.7 141.8 1 2
## 6 genuine 215.7 130.8 130.5 9.0 10.1 141.4 1 1
## hierarchical_agglomerative hierarchical_divisive gmm
## 1 1 1 1
## 2 2 2 2
## 3 2 2 2
## 4 2 2 2
## 5 2 2 2
## 6 1 2 2Nel blocco successivo, calcoliamo le metriche di valutazione come nei casi precedenti.
# K-means
ari_kmeans <- adjustedRandIndex(banknote$Status, banknote$kmeans)
mapping_kmeans <- bestMapping(banknote$Status, banknote$kmeans)
cm_kmeans <- mapping_kmeans$cm
acc_kmeans <- cm_kmeans$overall["Accuracy"]
sens_kmeans <- cm_kmeans$byClass["Sensitivity"]
spec_kmeans <- cm_kmeans$byClass["Specificity"]
distances_km <- t(apply(data_scaled, 1, function(x) {
apply(km_model$centers, 1, function(cent) sqrt(sum((x - cent)^2)))
}))
score_kmeans <- distances_km[,2] / rowSums(distances_km)
roc_obj_kmeans <- roc(banknote$Status, score_kmeans)
auc_kmeans <- auc(roc_obj_kmeans)
# K-medoids
ari_kmedoids <- adjustedRandIndex(banknote$Status, banknote$kmedoids)
mapping_kmedoids <- bestMapping(banknote$Status, banknote$kmedoids)
cm_kmedoids <- mapping_kmedoids$cm
acc_kmedoids <- cm_kmedoids$overall["Accuracy"]
sens_kmedoids <- cm_kmedoids$byClass["Sensitivity"]
spec_kmedoids <- cm_kmedoids$byClass["Specificity"]
medoids <- as.matrix(pam_model$medoids)
distances_kmed <- t(apply(data_scaled, 1, function(x) {
apply(medoids, 1, function(med) sqrt(sum((x - med)^2)))
}))
score_kmedoids <- distances_kmed[,2] / rowSums(distances_kmed)
roc_obj_kmedoids <- roc(banknote$Status, score_kmedoids)
auc_kmedoids <- auc(roc_obj_kmedoids)
# Gerarchico Agglomerativo
ari_hierAgg <- adjustedRandIndex(banknote$Status, banknote$hierarchical_agglomerative)
mapping_hierAgg <- bestMapping(banknote$Status, banknote$hierarchical_agglomerative)
cm_hierAgg <- mapping_hierAgg$cm
acc_hierAgg <- cm_hierAgg$overall["Accuracy"]
sens_hierAgg <- cm_hierAgg$byClass["Sensitivity"]
spec_hierAgg <- cm_hierAgg$byClass["Specificity"]
auc_hierAgg <- NA
# Gerarchico Divisivo
ari_hierDiv <- adjustedRandIndex(banknote$Status, banknote$hierarchical_divisive)
mapping_hierDiv <- bestMapping(banknote$Status, banknote$hierarchical_divisive)
cm_hierDiv <- mapping_hierDiv$cm
acc_hierDiv <- cm_hierDiv$overall["Accuracy"]
sens_hierDiv <- cm_hierDiv$byClass["Sensitivity"]
spec_hierDiv <- cm_hierDiv$byClass["Specificity"]
auc_hierDiv <- NA
# GMM
ari_gmm <- adjustedRandIndex(banknote$Status, banknote$gmm)
mapping_gmm <- bestMapping(banknote$Status, banknote$gmm)
cm_gmm <- mapping_gmm$cm
acc_gmm <- cm_gmm$overall["Accuracy"]
sens_gmm <- cm_gmm$byClass["Sensitivity"]
spec_gmm <- cm_gmm$byClass["Specificity"]
score_gmm <- gmm_model$z[,2]
roc_obj_gmm <- roc(banknote$Status, score_gmm)
auc_gmm <- auc(roc_obj_gmm)Infine, raccogliamo i risultati in un dataframe riassuntivo.
results_df <- data.frame(
Metodo = c("K-means", "K-medoids", "Gerarchico Agglomerativo", "Gerarchico Divisivo", "GMM"),
ARI = c(ari_kmeans, ari_kmedoids, ari_hierAgg, ari_hierDiv, ari_gmm),
Accuracy = c(acc_kmeans, acc_kmedoids, acc_hierAgg, acc_hierDiv, acc_gmm),
Sensitivity = c(sens_kmeans, sens_kmedoids, sens_hierAgg, sens_hierDiv, sens_gmm),
Specificity = c(spec_kmeans, spec_kmedoids, spec_hierAgg, spec_hierDiv, spec_gmm),
AUC = c(auc_kmeans, auc_kmedoids, auc_hierAgg, auc_hierDiv, auc_gmm)
)
print(results_df)## Metodo ARI Accuracy Sensitivity Specificity AUC
## 1 K-means 0.8830121 0.970 0.97 0.97 0.9987
## 2 K-medoids 0.8272389 0.955 1.00 0.91 0.9946
## 3 Gerarchico Agglomerativo 0.5160469 0.860 0.95 0.77 NA
## 4 Gerarchico Divisivo 0.9212040 0.980 0.98 0.98 NA
## 5 GMM 0.9212042 0.980 0.99 0.97 0.9983Rispetto ai risultati con il dataset completo o ridotto, l’uso di solo 3 variabili ottenute con PCA ha ridotto la performance globale, in particolare per K-medoids e Gerarchico Agglomerativo. Sebbene l’accuratezza resti elevata, l’ARI e la specificità mostrano una leggera diminuzione, il che suggerisce che la riduzione delle variabili ha compromesso la capacità di questi algoritmi di catturare la struttura dei dati. In particolare, l’algoritmo Gerarchico Agglomerativo mostra una performance molto più bassa in termini di ARI e accuratezza, indicando che l’algoritmo è più sensibile alla perdita di variabili discriminanti.
In conclusione, dal punto di vista della misurazione delle performance con tecniche esterne (avendo a disposizione la vera classe delle osservazioni), l’algoritmo GMM è risultato il migliore nei tre diversi dataset, anche se, tutto sommato, gli altri algoritmi non hanno ottenuto dei cattivi risultati.
7. Conclusioni
Questo progetto, realizzato a scopo didattico, ha analizzato diverse tecniche di clustering applicate al dataset banknote, un dataset semplice con due classi ben distinte. Inizialmente sono state utilizzate misure interne – quali il BIC, l’indice silhouette e altre metriche di compattezza – per valutare la qualità dei cluster senza ricorrere alle etichette reali. Successivamente, l’analisi è stata estesa con misure esterne (Adjusted Rand Index, accuratezza, sensibilità, specificità e AUC) per confrontare le assegnazioni dei cluster con la vera classe delle osservazioni.
I risultati mostrano che, quando si utilizzano tutte le variabili o solo quelle discriminanti (Bottom e Diagonal), algoritmi come GMM e, in alcuni casi, K-means e K-medoids raggiungono prestazioni eccellenti, con ARI e AUC quasi perfetti. Al contrario, i metodi gerarchici, in particolare quello agglomerativo, hanno ottenuto risultati leggermente inferiori, anche se in alcune configurazioni il clustering gerarchico divisivo si è dimostrato competitivo. L’approccio GMM risulta particolarmente robusto, grazie alla possibilità di estrarre score continui tramite le probabilità a posteriori.
Questi risultati, sia interni che esterni, evidenziano come una corretta selezione o trasformazione delle variabili (ad es. tramite PCA) possa semplificare il modello senza compromettere la capacità di distinguere correttamente le banconote genuine da quelle contraffatte.
Per il futuro, un possibile sviluppo potrebbe essere l’applicazione di questi metodi a dataset più complessi, nonché l’ottimizzazione degli score per metodi hard. Inoltre, integrare tecniche di soft clustering o combinare il clustering con approcci supervisionati potrebbe migliorare ulteriormente la generalizzazione e la robustezza dei modelli, offrendo una visione più completa delle performance in contesti reali.